# Generate Figure 5
y <- seq(1989,2022,by = 3)

# Calculate ILFs of results from "aggregate_inequality.R"
ilfs <- lapply(1:length(y),function(i){
  Reduce("+",lapply(1:5,function(k){agg[[i]][[k]]$estimate}))/5
})

# plot
sequence <- rep(seq(from = 0, to = 1, by = 1/(nrow(ilfs[[1]]) - 1)))
s <- 0
t <- 0
x <- rep(0,length(sequence)^2)
y <- x
for (i in 1:length(sequence)){
  s <- s+1
  k <- length(sequence)*(s - 1) + 1
  j <- length(sequence)*s
  x[k:j] <- rep(sequence[s],length(sequence))
  y[k:j] <- sequence
  
}

# Calculate the coordinates of the corners of the identical distribution
# There are more alphas than necessary.
alpha <- c(0.1,0.25,0.3,0.4,0.5,0.6,0.7,0.75,0.8,0.9,0.95,0.99)
egal <- rep(0,12)
for (k in 1:12){
  # This formula is given in paper section 1.2.2 example 1, set d = 2
  egal[k] <- uniroot(function(x){x*(1-log(x)) - alpha[k]},c(0.001,1))$root
}

lr <- data.frame(x,y,sapply(c(2,7,8,12),function(i){c(t(ilfs[[i]]))}))
colnames(lr) <- c("x","y",as.character(c(1992,2007,2010,2022)))
lr <- reshape2::melt(lr, id.var = c("x","y"))

levels <- c(0.6,0.8,0.95)
sz <- 1
sz2 <- 1.2

lcontours <- ggplot(data=lr,aes(x=x,y=y,z=value)) + 
  geom_contour(aes(colour = variable),breaks = levels, size =sz) + 
  annotate("text", x=0.63, y=0, label= "Income", family = 'serif',size = 5.5) +
  annotate("text", x=0, y=0.63, label= "Wealth", family = 'serif',size = 5.5,angle = 270) +
  annotate("label", x=0.9, y=0.58, label= "0.95", family = 'serif',size = 5.5) +
  annotate("label", x=0.9, y=0.40, label= "0.90", family = 'serif',size = 5.5) +
  annotate("label", x=0.9, y=0.18, label= "0.80", family = 'serif',size = 5.5) +
  annotate("label", x=0.9, y=0.07, label= "0.60", family = 'serif',size = 5.5) +
  annotate("segment",x = 0, y = 0, xend = 0.63, yend = 0.63,linetype = 3, size=sz)+
  scale_color_manual(values=c("1992" = "black","2007"="#ddaa33","2010"="#994455","2022" = "#004488")) +
  coord_cartesian(xlim = c(0, 0.65),ylim = c(0, 0.65)) +
  theme(plot.title = element_text(hjust = 0.5)) + xlab("")+ ylab("")+ 
  guides(col= guide_legend(title= "Year"))+
  theme(panel.border = element_rect(colour = "black", fill=NA, size=0.6),
        panel.background = element_rect(fill = 'transparent'),
        panel.grid = element_blank(),
        legend.position = c(0.85,0.85),
        legend.key = element_rect(fill = 'transparent'),
        text = element_text(family = "serif",size=22))+
  scale_y_continuous(breaks = seq(0,1,0.2),position = "right",labels = seq(0,1,0.2)) +
  scale_x_continuous(breaks = seq(0,1,0.2),position = "top",labels = seq(0,1,0.2)) 
ggsave(plot = lcontours, width = 6, height = 6, filename = "fig5.png")
